Quantum path-integral study of the phase diagram and isotope effects of neon 
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The phase diagram of natural neon has been calculated for temperatures in the range 17-50 K 
and pressures between 10~^ and 2 x 10'' bar. The phase coexistence between solid, liquid, and 
gas phases has been determined by the calculation of the separate free energy of each phase as a 
function of temperature. Thus, for a given pressure, the coexistence temperature was obtained by 
the condition of equal free energy of coexisting phases. The free energy was calculated by using 
non-equilibrium techniques such as adiabatic switching and reversible scaling. The phase diagram 
obtained by classical Monte Carlo simulations has been compared to that obtained by quantum 
path-integral simulations. Quantum effects related to the finite mass of neon cause that coexistence 
lines are shifted towards lower temperatures when compared to the classical limit. The shift found 
in the triple point amounts to 1.5 K, i.e., about 6% of the triple-point temperature. The triple-point 
isotope effect has been determined for ^"Ne, ^^Ne, '^^Ne, and natural neon. The simulation data 
show satisfactory agreement to previous experimental results, that report a shift of about 0.15 K 
between triple-point temperatures of ^°Ne and ^^Ne. The vapor-pressure isotope effect has been 
calculated for both solid and liquid phases at triple-point conditions. The quantum simulations 
predict that this isotope effect is larger in the solid than in the liquid phase, and the calculated 
values show nearly quantitative agreement to available experimental data. 
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I. INTRODUCTION 

The electronic structure of atoms and molecules, en- 
countered either in gas or condensed phases, is, within 
the framework established by the Born-Oppenheimer 
(BO) approximation, independent of the nuclear mass. 
Thus the BO potential energy surface (PES), that repre- 
sents the energy landscape for the motion of the nuclei, 
results to be independent of the isotope mass. Many 
processes and properties depend on the movement of the 
atoms on the PES. Then, isotope effects may exist be- 
cause the kinetic energy of the nuclei is mass dependent 
according to the laws of quantum mechanics. Thus, al- 
though potential energy parameters are isotope indepen- 
dent, the kinetic energy parameters are not.^ Phase co- 
existence is an equilibrium property that depends on the 
kinetic energy of the atoms (think of a liquid in equi- 
librium with its vapor). Hence, the phase diagram of a 
given substance should depend on its nuclear mass and 
then display an isotope effect. We stress that this be- 
havior is strictly of quantum mechanical nature, i.e., the 
classical limit in statistical thermodynamics incorrectly 
predicts that phase coexistence is a mass independent 
property. This erroneous prediction can be considered as 
a consequence of the equipartition theorem for the kinetic 
energy, that establishes that in the classical limit this en- 
ergy is only a function of temperature and is therefore 
independent of the atomic mass. 

The properties of rare-gas solids or fluids usually 
change as a function of the atomic mass in a systematic 
fashion, that has been used to study effects of the zero- 
point energy and isotope mass. There have been many 
studies of quantum effects in the phase coexistence of 
neon, which were initially motivated by the observation 



that neon deviates from the principle of corresponding 
states, that is however rather accurately followed by ar- 
gon, krypton, and xenon,^ The calculation of quantum 
corrections to the classical limit of the phase diagram 
of neon is thus interesting because allows us to quantify 
systematic errors of treating the neon atoms as particles 
moving classically in the PES. However, we stress that 
these quantum corrections do not represent any kind of 
measurable property, i.e., they can not be directly com- 
pared to experimental data, because there is no way to 
perform measurements of the phase behavior in the clas- 
sical limit. In this respect, the calculation of isotope 
effects in the phase diagram of neon results specially in- 
teresting because they can be directly compared to exper- 
imental data, offering then a direct test of the capability 
of the theoretical approach. Neon has three stable iso- 
topes (^"Ne, ^^Ne, and ^^Ne) and the experimental study 
of isotope effects in its phase diagram is historically in- 
teresting. The observation by Keesom and Dijk^ in 1931 
that ^"Ne was more volatile than ^^Ne at temperatures 
just above the triple point was the first demonstration of 
an isotope effect on the physical or chemical properties of 
matter, and this experiment was considered to be a proof 
of the existence of the zero-point energy^i An accurate de- 
termination of the vapor pressure difference between nat- 
ural neon and pure ^°Ne and ^^Ne isotopes has been made 
by Bigeleisen and Rotbi and Furokawaj^ The vapor pres- 
sure isotope effect (VPIE) amounts to about 5% of the 
equilibrium pressure, being slightly larger for solid than 
for liquid phases. The triple-point isotope effect (TPIE) 
has been also measured for the ^"Ne and ^^Ne isotopes?^ 
and amounts to a temperature shift of 0.15 K. This differ- 
ence might seem small, but it is relevant for thermometry 
as the triple point of natural neon (Tp — 24.553 K) is used 
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as fixed-point temperature standard in the International 
Temperature Scale of 1990 (ITS-90), and any deviation 
of the isotopic composition of natural neon leads to am- 
biguity in the reference temperatureii 

The most employed quantum approaches to the phase 
diagram of neon are based either on perturbational meth- 
ods using the Wigner-Kirkwood expansion of the quan- 
tum partition function in powers of h, or on computer 
simulations using the path integral (PI) formulation of 
statistical mechanics. We are not aware of PI simula- 
tions of isotopic effects in the phase coexistence of neon. 
However, PI MC simulations of isotope effects in the lat- 
tice parameter and heat capacity of solid neon have been 
reported)^ showing good agreement to experimental data. 
The VPIE of liquid neon has been recently studied^ using 
a classical integral equation theory in order to obtain the 
parameters required for the statistical theory of isotope 
effects developed by Bigeleisen in 196lJ^ The isotopic 
shift in the helium melting curve has been studied by PI 
MC simulations at temperatures where nuclear exchange 
effects associated to the nuclear quantum statistics are 
negligible. ^^'^^ Isotope effects in the solid-liquid coex- 
istence of Lennard- Jones (LJ) solids have been studied 
by the method of generating Landau free energy curves. 
This method was applied to study the melting curves 
of solid H2 and D2, for which quantum effects are ex- 
pected to be much larger than for neon, due to its lighter 
molecular mass. At about 12 MPa the calculated melting 
temperature of ortho-Y)2 and para-Yi.2 was of 22.3 and 18 
K, respectively.^'^ 

In a summary of theoretical calculations of quantum 
corrections in neon, we should mention the study by 
Hansen and Weis^ of the coexistence curve of neon by 
a Wigner-Kirkwood expansion to order fi^. They found 
a lowering of the triple-point temperature of about 1.4 
K and of the critical temperature of about 1.8 K due to 
quantum corrections. The quantum effect at the liquid- 
gas critical point was studied by Young by the Feynman- 
Hibbs effective potential method. An interesting pre- 
diction was that quantum effects modify the critical pres- 
sure, volume, and temperature, but no effect was pre- 
dicted on the critical exponents. The neon liquid-vapor 
interface at the experimental triple point was studied by 
PI Monte Carlo (MC) simulations and compared to clas- 
sical simulations, showing that the quantum interface is 
approximately one quarter of interatomic spacing wider 
than the classical counterpart fiS The liquid- vapor equi- 
librium of neon has been studied by Gibbs ensemble MC 
calculations by both the PI formulation,^^ and a Wigner- 
Kirkwood expansionJ^^ Quantum effects in the liquid- 
vapor equilibrium were shown to be important and more 
significant in the liquid state. The critical temperature in 
the quantum simulation was found about 8% lower than 
in the classical limit This value results to be larger 
than that of 4% predicted by a Wigner-Kirkwood expan- 
sion in Ref. [lj|. The solid-liquid coexistence line of 
neon has been studied using a Wigner-Kirkwood effec- 
tive potential at pressures between 10^ and 4 x 10'^ MPa, 



showing that quantum effects are negligible in this range 
of pressures A recent PI MC simulation has shown that 
at ambient pressure the melting temperature of neon is 
lowered by about 6 % (i.e. 1.5 K) when quantum effects 
are included 1^ 



In the present work we study the phase diagram of 
natural neon for temperatures in the range 17-50 K and 
pressures between 10~^ and 2 x 10'^ bar by calculation 
of the phase coexistence between solid, liquid, and gas 
phases. Moreover the VPIE of solid and liquid phases at 
triple-point conditions and the TPIE of ^"Ne and ^^Ne 
are also calculated and compared to experimental data. 
All results are obtained by PI MC simulations. Phase co- 
existence is obtained by requiring equality of the Gibbs 
free energy, G, of the phases in equilibrium. The calcu- 
lation of G is performed using the Adiabatic Switching 
(AS)^^ and Reversible Scaling (RS)^ approaches that are 
based on algorithms where the Hamiltonian (or a state 
variable as the pressure or inverse temperature) changes 
along the simulation run. The capability of both AS and 
RS methods to calculate free energies in the context of 
PI simulations has been recently analyzed in the study of 
neon melting, and we will make extensive use of these 
free energy techniques in the present investigation. The 
influence of quantum effects in the phase diagram of neon 
will be quantified by comparing the predicted phase di- 
agram calculated with and without atomic quantum ef- 
fects included. The interatomic interaction between Ne 
atoms will be described by a Lennard- Jones pair po- 
tential, in analogy to our previous studies We 
recognize that a better option for the neon pair poten- 
tial is given by the fit to ab initio electronic structure 
calculations r^iS^i^ however we have chosen the simple 
LJ pair potential for several reasons: the phase diagram 
in the classical case can be checked against a large num- 
ber of previous LJ simulations,—"— quantum effects in 
the three coexistence lines (solid-gas, solid-liquid, and 
liquid-gas) has not yet been calculated by PI simulations 
using the same LJ potential parameters, and the isotope 
effects in the phase diagram depend more on the kinetic 
energy than on the potential energy parameters. 



The structure of this paper is as follows. In Sec. Illlwe 
present the computational conditions employed in the PI 
simulations and the technique used to evaluate the free 
energy as a function of the isotope mass. The results 
obtained by PI and classical simulations for the three 
coexistence lines of neon are presented in Sec. IIIII The 
phase diagram of neon is displayed both in the pressure- 
volume (P — V) and in the density-temperature (p — T) 
domains. The TPIE of 2°Ne and ^^Ne is calculated and 
compared to experimental data in Sec. IIVI while the 
VPIE of both solid and liquid phases of neon are derived 
in Sec. [Vj including also the comparison to experiment. 
Finally, we summarize our conclusions in Sec. IVIl 
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II. METHODOLOGY 
A. Computational conditions 

In the PI formulation of statistical mechanics the par- 
tition function is calculated through a discretization of 
the integral representing the density matrix. This dis- 
cretization serves to define cyclic paths composed by a 
finite number L of steps, that in the numerical simu- 
lation translates into the appearance of L replicas (or 
beads) of each quantum particle. Then, the implementa- 
tion of PI simulations relies on an isomorphism between 
the quantum system and a classical one, derived from 
the former by replacing each quantum particle (here, 
atomic nucleus) by a ring polymer of L classical parti- 
cles, connected by harmonic springs with a temperature- 
and mass-dependent force constant. Details on this com- 
putational method are given elsewhere<^"— In this work 
the technical implementation of the PI simulations in the 
NPT and NVT ensembles {N being the number of par- 
ticles) is identical to that presented in Ref. (20| . 

For the PI simulation of the phase diagram of natu- 
ral neon, the rare gas atoms were treated as particles of 
mass mo = 20.18 amu interacting through a LJ poten- 
tial with parameters e = 3.2135 meV and a = 2.782 A. 
The value of the parameter e differs from that employed 
in our study of the melting transition of neon at ambi- 
ent pressure, where we found that the classical melting 
temperature {T„i^ci =24.95 K) was closer to the exper- 
imental data (r„i =24.55 K) than the quantum result 
{Tm.q =23.43 K). Note that the classical melting tem- 
perature, Tm.ci, is proportional to the value of e, and a 
significant deviation of Tm^q from this proportionality is 
not expected in the quantum case. Thus, the parameter 
e used in Ref. 0] has been increased by a scaling factor 
of a =1.042, which is roughly the factor needed to match 
experimental and calculated T^.q temperatures. The LJ 
potential was truncated for distances larger than Tc =4f7 
and shifted to zero to perform the MC sampling. The 
internal energy was corrected by adding back this poten- 
tial shift in order to approximate the properties of the 
full LJ potential. 

PI MC simulations were performed on cubic cells con- 
taining 500 atoms, assuming periodic boundary condi- 
tions. For the liquid, standard long-range corrections 
were computed assuming that the pair correlation func- 
tion is unity, g{r) = 1 for r > Tc, leading to well-known 
corrections for the pressure and internal energies <^ In the 
solid, the long range correction was a static-lattice ap- 
proximation (called LR-s in Ref. [1^), that implies that 
atoms beyond a certain number of neighbor shells are ar- 
ranged as in the perfect crystal ordering. From the value 
of the cutoff radius, we determine that the first twelve 
neighbor shells (representing 248 neighbors of a central 
atom) were dynamically included in the simulation. This 
approximation leads also to simple long-range corrections 
for the pressure and internal energiesj^i^ The virial esti- 
mator was used for the calculation of the kinetic energyj^ 



Solid phase simulations were started with the atoms on 
their ideal face-centered-cubic (fee) positions, while the 
liquid and gas phase simulations were initialized by us- 
ing random positions. The number of beads was set as 
L = 15, which has been shown to be adequate for the 
temperatures and pressures employed in this studyi^ '^'^i^^ 
Typical runs consisted of 3 x 10"* MC steps (MCS) for 
equilibration, followed by at least 10^ MCS for the calcu- 
lation of equilibrium properties. Some simulations were 
performed with larger runs up to 6 x 10^ MCS, to check 
the convergence of the results. Classical results were ob- 
tained by setting the number of beads L—1 in our PI MC 
codes. 



B. Free energy calculations 

The thermodynamic integration (TI) is an standard 
technique that allows us to obtain the free energy of a 
system by the calculation of the reversible work needed 
to change the original system into a reference state of 
known free energy.—"— The AS method is an alternative 
to the standard TI, where the reversible work is obtained 
by slowly changing the system Hamiltonian along the 
simulation run.^^ RS is an efficient technique to obtain 
the free energy as a function of a state variable, typi- 
cally T or PM- The reversible work is calculated here, in 
way similar to the AS method, by slowly changing the 
state variable along the simulation run. The PI imple- 
mentation of both AS and RS methods has been recently 
described in the study of solid-liquid coexistence of neon 
at ambient pressure and will not be repeated herejSS An 
additional free energy method employed in this work is 
the Morales and Singer (MS) approach, that allows us 
to obtain the free energy difference between a quantum 
system and its corresponding classical limit in either the 
NVT or NPT ensembles,^!^ Here we merely summa- 
rize how these methods fit into our calculation procedure 
of the phase coexistence of neon. 

The first step is to obtain the free energy of each phase 
(solid, liquid, or gas) at some appropriate reference point 
in the (T, P) plane. An NPT simulation of each phase 
at the chosen reference point is run to determine its equi- 
librium volume. In the case of the solid, the Helmholtz 
free energy, F, is then calculated by an AS simulation in 
the NVT ensemble, using the Einstein crystal as a refer- 
ence. For the quantum case the reference is a quantum 
Einstein crystal^^ while in the classical limit a classical 
Einstein crystal is used as reference.'*^ 

In the case of a quantum fluid (either liquid or gas) , 
the free-energy difference between the quantum and clas- 
sical fluid is calculated with the MS approach in either 
the NVT or NPT ensembles.^O.^i The excess free energy 
of the classical fluid is then obtained from the Johnson 
et al. parameterizationi^ Finally, the ideal gas contri- 
bution is added to obtain the total free energy of the 
quantum fluid 

In the case of the solid and liquid phases a finite size 
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correction was applied to the free energy Fn obtained for 
N —500 atoms by using the foUowing relation 



F = 



N 



a{T) 
i 



(1) 



where F represents the infinite-size free energy, a(T) is a 
constant that depends linearly on temperature (see Ref. 
[23]) and the scaling factor a — 1.042 takes into account 
the change in the e parameter employed in this work and 
in Ref. [20|. The value of a(T) for both solid and liquid 
phases of neon is obtained from Tab. II and Eq. (21) in 
Ref. 

The coexistence temperature between two phases is 
then calculated at a given pressure by a free energy ap- 
proach, by which the free energy of both phases is calcu- 
lated over an overlapping range of temperatures by the 
RS method. The coexistence temperature is then deter- 
mined by the point where the Gibbs free energy of both 
phases is identical. 



C. Free energy isotope effect 



where Am = (m — m[))/{M — 1). The free energy differ- 
ence, AF = F{m) — F{mo), is then obtained as 



M 



af^-j: 



K^{mi) 
Wi Am 



(6) 
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where Ky^rrii) is the virial estimator of the kinetic energy 
at step i. Wi is an integration weight factor (wi — 1 ex- 
cept for i=l or i = N, where Wi =0.5). A calculation of 
isotope effects on solubilities of hydrogen has been done 
by a Fourier PI relation equivalent to Eq. ([4]) by using 
TL— In the case of the NPT ensemble the mass depen- 
dence of the Gibbs free energy can be derived from the 
following relation analogous to Eq. 



G(m) = G(mo) 



Kim'] 



dm' . 



(7) 



For the calculation of VP IE at a given temperature it is 
useful to determine the Gibbs free energy as a function of 
pressure, provided that an initial value, G{Po), is already 
known. We have employed the following equation 



The isotope effect in the phase coexistence of a sub- 
stance is determined by the free energy difference be- 
tween isotopes. Given an isotope of mass m, and assum- 
ing that the Helmholtz free energy of the natural isotopic 
composition, i^(mo), is known by any of the methods 
summarized in Sec. Ill Bl then the unknown free energy, 
F{m), can be calculated by the expression 



m 

F(m) = F(mo) + J 



^EMldm' 
dm' 



(2) 



mo 



By considering the following thermodynamic definition 
of the kinetic energy K , 



we get 



K{m') = —m' 



F{m) = F(mo) 



, dF{m') 
dm' ' 



K{m') 



m' 



dm' 



(3) 



(4) 



I' 

G(P) = G(Po) + J VdP' 



(8) 



where we have made used of the relation V = dG/dP. 
The calculation of the integral in the last equation is also 
done by the RS method, i.e., the pressure P' is uniformly 
changed between Pq and P along the simulation run and 
the value of G(P) is calculated "on the fly" by solving 
numerically the integral in Eq. ([8]) as the simulation 
proceeds 1^ 



III. PHASE DIAGRAM OF NEON 

The phase coexistence of natural neon has been de- 
rived from PI MG simulations. We present the results ob- 
tained for the three coexistence lines (solid-liquid, solid- 
gas, and liquid-gas). The studied temperature range is 
between 17-50 K and pressures were varied between 10~^ 
and 2 x 10'^ bar. The quantum results are compared to 
the classical limit and to available experimental data. 



This expression shows that the kinetic energy is the mag- 
nitude that determines the free energy difference between 
two isotopes. The implementation of the last equation in 
our PI simulations has been done by the AS method. 
Thus, in a simulation run with a number M of MGS, the 
particle mass is changed uniformly from the initial value 
mo to the final value m. The actual value of the mass at 
the simulation step i is. 



mi — mo + (i — l)Am 



(5) 



A. Solid-liquid coexistence 

The free energy of both solid and liquid phases of neon 
over an overlapping range of temperatures that encom- 
passes the melting point is shown in Fig.[l]for P— 1 bar. 
The classical melting point is found at 26.1 K. The main 
quantum effect in the free energy is a shift of the classical 
data towards higher values. As this energy shift results 
larger for the solid than for the liquid phase, the melting 
temperature is displaced downward by about 1.5 K. The 




Figure 1: Melting temperature for fee neon at 1 bar obtained 
by the condition of equal Gibbs free energy of the coexisting 
phases. Results are shown for both quantum and classical 
descriptions of solid (s) and liquid (1) phases. 



slopes of the free energy curves are proportional to -S {S 
being the entropy). Prom the figure, it is obvious that the 
entropy of the liquid phase is larger than that of the solid. 
Less obvious to the eye is that the entropy of both phases 
increases when quantum effects are included. At melting 
conditions (P=l bar, Tjn,ci = 26.1 K, and T^,, = 24.6 
K) the entropy of the solid phase amounts to 14.7 J/mol 
K (quantum) and 10.2 J/mol K (classical), while for the 
liquid phase the corresponding values are 27.7 J/mol K 
(quantum) and 23.8 J/mol K (classical). Therefore the 
consideration of quantum effects implies an increase of 
44 % in the solid entropy and of 16 % in the liquid one. 

The solid-liquid phase consistence has been calculated 
for 9 different pressures up to 2 Kbar and the result for 
the melting line is shown in Fig. 2. We observe that the 
quantum melting line is displaced with respect to the 
classical limit. At zero pressure the quantum shift in the 
melting temperature, T^, amounts to 1.6 K (about 6 % 
of Tm), while at 2 Kbar this shift is only of about 0.6 
K (about 1 % of Tjn)- Therefore, the quantum effect in 
the melting temperature decreases as pressure goes up. 
This fact is in line with the Wigner-Kirkwood analysis 
of Ref. [l^ for pressures between 10'^ and 2 x lO''^ bar. 
Our PI MC data show a reasonable overall agreement to 
experimental results^i The agreement is in part origi- 
nated by our selection of the LJ parameter e to fit the 
experimental melting temperature (24.55 K) at ambient 
pressure (see Subsec. Ill A[) . The broken line in the figure 
shows the classical simulation results of Ref. 27] for the 
LJ parameters employed here. These data are slightly 
shifted with respect to our classical melting line. We ob- 
tain somewhat higher melting temperatures for a given 
pressure. A possible reason for this discrepancy is that in 
both simulations the long range correction used for the 



Figure 2: Solid-liquid (s-1) coexistence line of neon in the 
P — T domain. The continuous line is the experimental result 
from Ref. [3 ]• Open circles and squares are the results of 
our PI MC and classical MC simulations, respectively. The 
broken line is the empirical fit to classical simulations given 
in Ref. 

solid phase and the LJ potential cutoff arc different. 



B. Solid-gas coexistence 

The Gibbs free energy G of both solid and gas phases 
of neon is shown as a function of temperature in Pig. |3l 
The results derived for the gas from the Johnson parame- 
terization are shown for five pressures between 10~^ and 
0.5 bar. Neither quantum effects nor deviations from 
ideal behavior are significant in the displayed range of 
pressures for the gas free energies. The free energy of the 
solid corresponds to the data calculated at a pressure of 
1 bar. At the scale of this figure, solid free energies at 
lower pressures are indistinguishable from the P = 1 bar 
data and are not represented (the PV term of the Gibbs 
free energy is of the order of 1 J/mol at these pressures). 
However, lower pressure results were calculated for the 
solid phase to determine the solid-gas coexistence line. 
In contrast to the gas behavior, quantum effects are im- 
portant in the solid, and thus both classical and quantum 
free energy results are plotted in Pig. [3l It is interesting 
to note that at temperatures about 19 K the classical 
free energy of the solid displays a maximum, therefore 
below this temperature the entropy of the classical solid 
becomes negative and in fact it would tend to — cx) as tem- 
perature goes down to zero. A simple way to derive this 
result is to work out analytically a harmonic model for 
the classical solid. This behavior of the classical entropy 
is known to be incompatible with the Nerst theorem that 
constitutes the third principle of thermodynamics. 

The crossing points of the solid and gas free energy 
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Figure 3: Gibbs free energies of the solid (s) and gas (g) 
phases of neon as a function of temperature. Gas results are 
shown for five different pressures. The free energy of the solid 
corresponds to a pressure of 1 bar. For the solid, the value of 
G for pressures between and 1 bar are indistinguishable at 
the scale of the figure. Results for both quantum and classical 
simulations of solid fee neon are shown. 
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Figure 4: Solid-gas (s-g) coexistence curve of neon in the P—T 
domain. The continuous line is the experimental result from 
Ref. Open circles and squares are the results of our PI 

MC and classical MC simulations, respectively. The dotted 
line is a polynomial fit of fourth degree to the PI MG results. 
The broken line is a classical fit obtained by an empirical EoS 
in Ref. [H. 



curves, G{T), at a given value of P determine the subli- 
mation curve of the sohd. We have calculated the solid- 
gas coexistence line at eighth values of P in the range be- 
tween 10~^ and 0.8 bar. The results are shown in Fig. |3] 
and compared to available experimental data up to the 
triple-point temperature of neoniiS The PI MC results 
(open circles) display a constant shift of about 0.6 K 
with respect to the experimental data. We note that by 
the same scahng considerations as presented in Subsec. 
Ill Al for the LJ parameter e, a decrease in the value of 
e would produce a shift of our calculated solid-gas coex- 
istence Hne toward lower temperatures. In this way we 
might improve the agreement of our sublimation line to 
experiment. However, the calculated melting line would 
then appear also shifted toward lower temperatures, and 
thus the good agreement to experimental data shown in 
Fig. [2] would be lost. We conclude that a simple LJ pair 
potential seems to be unable to describe both solid-gas 
and solid-liquid coexistence lines of neon in good agree- 
ment to experiment. 

Our classical result for the sublimation line of neon is 
also given in Fig. The shift with respect to the PI 
MC results is entirely due to quantum effects of the solid 
phase, as the gas behaves as an ideal gas for the displayed 
values of pressure and temperature. The broken line in 
the figure shows the fit of van der Hoef of the solid-gas 
coexistence using an empirical equation of state (EoS) 
for the LJ solid and the Johnson parameterization for 
the gas4S. Our classical results are consistent with this 
fit. 
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Figure 5: Determination of the boiling point, Tt, of neon at 
20 bar. The calculated G{T) curves for liquid (1) and gas (g) 
phases are displayed as a function of temperature. Results 
are displayed for both quantum and classical simulations. 



C. Liquid-gas coexistence 

The Gibbs free energy curves, G{T), of liquid and gas 
neon phases at pressure of 20 bar are displayed in Fig. [S] 
At this pressure quantum corrections to the free energy of 
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Figure 6: Liquid-gas (1-g) coexistence curve of natural neon 
in the P — T domain. The continuous line is the experimental 
result from Ref. [46|. Open circles and squares are the results 
of our PI MC and classical MC simulations, respectively. The 
dotted line is a polynomial fit of sixth degree to the PI MC 
results. The broken line is the coexistence line obtained by 
using the parameterized EoS of Johnson.-2i 



the gas are significant, although smaller than in the case 
of the liquid phase. The shift in the boiling temperature 
due to quantum effects amounts to 1.4 K. The classical 
G(T) curve for the gas phase was derived from the John- 
son parameterization,'^^ while for the liquid phase the 
classical curve was obtained by the RS method started 
at a temperature of T = 44 K. We have proceeded in this 
way to check the consistency of our classical MC liquid 
simulations against the Johnson EoS."^^ 

The liquid-gas coexistence line of natural neon is dis- 
played in Fig. ini Our classical results, derived from RS 
simulations of the free energy in the liquid phase, are con- 
sistent with the coexistence line directly obtained from 
the Johnson EoS,^ The PI MC boiling line is shifted 
with respect to the classical result. Quantum corrections 
in the gas phase were found to be significant at pres- 
sures above 5 bar. Our quantum results deviate from 
the experimental boiling line4^ Similar differences have 
been reported in Gibbs ensemble MC simulations of the 
vapor-liquid phase coexistence of neoufii and the main 
reason for this discrepancy seems to be related to the ne- 
cessity of including three body terms in the interaction 
potential 

We have performed NPT simulations of both liquid 
and vapor phases of neon at coexistence conditions to ob- 
tain the saturation densities. The critical temperature, 
Tc, has been then estimated by fitting the calculated den- 
sities of the liquid and vapor phases, pi and pg, to the 
scaling law^ 
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density (gr/cm ) 



Figure 7: Density of coexistence of vapor and liquid phases. 
Experimental data are taken from Ref. [i^ ], and classical 
results from Ref. [sj ]. The broken line is a fit of our PI MC 
results (open circles) to the scaling law presented in the text 
[Eqs. (|9} and (I10|l ]. Closed circles correspond to the critical 
point in the different approaches. 



Table I: Comparison of the calculated critical point of neon 
with experimental data. 



experiment" PI MC classical' 



Tc{K) 
Pc(bar) 
Pc(g/cm^) 



44.4 
26.5 
0.483 



46.7 
27.7 
0.477 



49.0 
30.9 
0.483 



"From Ref |46l]. 
''Johnson EoSJi 



where = 0.32 is the non-classical critical exponent, 
which is expected to be valid also in the quantum case,— 
and 6 is a constant. The values of Tc and b where obtained 
by a least-squares fit of our PI MC results. The critical 
density, p^-, was estimated by a fit to the law of rectilinear 
diameters 4^ 



1 



{p,+Pg)=Pc + A{T-Tc) 



(10) 



Pi- Pg^ b{Tc - T) 



(9) 



where the constant A and the density, pc, are the fitted 
parameters. The PI MC estimation of the critical pres- 
sure, Pc, was obtained from the fit of the boiling line in 
Fig. [6] at temperature Tc- The critical parameters esti- 
mated from our quantum simulations are compared to 
the experimental data and to the Johnson EoS in Tab. U 
and the results for the saturation densities are presented 
in Fig. [71 The quantum correction to the saturation den- 
sity is of opposite sign for vapor and liquid phases, while 
the vapor coexistence density increases, the liquid one 
decreases. Although the PI MC saturation densities are 
closer to experiment than the classical result, the devia- 
tion from experiment is notorious, especially as one ap- 
proaches the critical point. This discrepancy is also found 
in previous simulations of the vapor-liquid coexistence of 
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Figure 8: Phase diagram of natural neon in the pressure- 
temperature domain. Results are shown for the three coex- 
istence lines between solid (s), liquid (1), and gas (g) phases. 
Continuous lines are derived from the PI MC results shown 
in Figs. [21 |4l and [S] while dotted lines are derived from the 
corresponding classical data in the same figures. Open circles 
represent the triple point and closed ones the critical point. 



Table III: PI MC results for thermodynamic properties of nat- 
ural neon at triple-point conditions. 

G (kJ/mol) PV (kJ/mol) U (kj/mol) S (J/mol K) 

solid -2.0492 0.0004 -1.6883 14.68 

liquid -2.0492 0.0005 -1.3719 27.57 

gas -2.0492 0.2002 0.2989 103.80 



For both quantum and classical results we find that the 
three coexistence lines meet within a small temperature 
interval of 0.01 K and within a pressure interval of 0.002 
bar. This small dispersion in the triple point determi- 
nation is an internal check of the consistency of our re- 
sults. Calculated triple-point properties are summarized 
in Tab. HI The value of Tjp predicted by the PI MC 
results is in good agreement to experiment, while classi- 
cal simulations result in a higher value of Ttp by about 
1.5 K. The equilibrium densities of the solid and liquid 
phases calculated by PI MC are also in good agreement 
to experiment. The calculated triple-point pressure, Ftp-, 
is lower than the experimental data, and we do not ob- 
serve any appreciable quantum correction to the value of 
Ptp. The calculated values of G, PV, internal energy (J7), 
and S at triple-point conditions are summarized in Tab. 

nil 



Table II: Comparison of predicted triple-point location of nat- 
ural neon with experimental data. 

experiment PI MC classical 



Ttv (K) 


24.553" 


24.55 


26.04 


Ptp (bar) 


0.433" 


0.315 


0.323 


ptp.s (g/cm^) 


1.424*' 


1.429 


1.494 


Ptp.i (g/cm^) 


1.247'' 


1.254 


1.313 


ptp,a(g/cm^) 




0.0032 


0.0031 



"Ref. 
''Ref. 



neon using different pair potentials and quantum correc- 
tions, and reflects the importance of three body terms in 
the interaction potential to correctly describe phase co- 
existence near critical point conditions! ^''''^^i^^ From the 
difference between calculated and experimental values of 
Tc in Tab. U one can infer that near critical point con- 
ditions both quantum corrections and three body terms 
should be of similar importance. 



E. T-p domain . 

NPT simulations of coexisting phases have been per- 
formed along the three coexistence lines determined in 
the previous Subsections to obtain the phase diagram of 
natural neon in the temperature-density domain. The 
PI MC and classical results are presented in Fig. [HI The 
largest deviation between classical and quantum results 
is found for the solid phase, especially at temperatures 
below the triple point. We observe that quantum correc- 
tions to the liquid coexistence density at temperatures 
between the triple (Ttp) and critical points (Tc) results 
to be larger for the liquid-gas line than for the liquid- 
solid one. This behavior is probably related to the trend 
that quantum effects are comparatively less important as 
pressure increases,-^ The liquid at temperatures between 
Ttp and Tc is under higher pressure if it is in equilibrium 
with the solid than if it is in equilibrium with its vapor 
(see the P — T phase diagram in Fig. [S] and note that 
the coexistence pressure increases in both cases due to 
quantum effects). 



D. P-T domain 

The three calculated coexistence lines of natural neon 
are presented in Fig. |5]in the P — T domain. Note that 
the pressure is plotted using a logarithmic scale. Contin- 
uous lines were derived from PI MC simulations, while 
dotted lines are classical results. The triple point is de- 
termined by the crossing of the three coexistence lines. 



IV. TRIPLE-POINT ISOTOPE EFFECT 

We turn now to the study of the TPIE. We have pre- 
viously obtained that the triple-point pressure, Ptp, of 
neon is nearly the same in both quantum and classical 
limits. Therefore one can reasonably assume that Ptp is 
nearly independent of the isotope mass. This fact sim- 




density (g/cm ) 



Figure 9: Phase diagram of natural neon in the temperature- 
density domain. Closed squares and circles are results ob- 
tained from classical MC and quantum PI MC simulations, 
respectively. Classical vapor and liquid densities were derived 
from the Johnson EoS,''* while the continuous line through 
the quantum vapor and liquid densities is the fit to the scal- 
ing law presented in the text [Eqs. @ and pop ]. Results are 
shown for the three coexistence lines between solid (s), liquid 
(1), and gas (g) phases. 



Table IV: Free energy differences (AGm, in J/mol) between 
neon isotopes and natural neon calculated at the triple point 
of natural neon. The calculated and experimental triple-point 
temperatures for each isotope are collected in the last lines. 

"°Ne ^^Ne ^^Ne 
AG,„ (s) 3.8 -16.9 -36.4 
AG™ (1) 3.6 -16.1 -34.8 

Ttp(K) (PI MC) 24.538 24.608 24.673 
Tfp(K) (exper.)° 24.540 - 24.687 

"Ref. [0. 

plifies the calculation of the TPIE that will be derived 
from a phase coexistence calculation at fixed pressure 
{Ftp = 0.315 bar). We have arbitrarily chosen to cal- 
culate the coexistence between solid and liquid phases, 
although any other two phases would be equally appro- 
priate. The melting point of the neon isotopes will be 
determined by the same procedure shown in Fig. [1] for 
natural neon. 

We take {Ttp,Ptp) of natural neon as reference point 
to calculate the Gibbs free energy of the solid and liquid 
phases of the neon isotope with mass m. The Gibbs free 
energy difference, 

AG„, = G(m) - G(mo) , (11) 

between the isotope of mass m and natural neon (toq) 
is derived by an AS calculation using Eq. ([7]). The free 
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Figure 10: Triple-point temperature differences for the differ- 
ent neon isotopes as a function of the variable x defined in the 
text [Eq. (|12[) ]. PI MC results are compared to experimental 
data. 

energy G(mo) was given in Tab. IIIIl and the results of 
AGm for both solid and liquid phases are summarized 
in Tab. IIVI as a function of the isotope mass m. The 
triple-point temperatures, T^p, obtained as the melting 
point of the isotopes at Ftp = 0.315 bar are also given 
in Tab. IIVI The difference between the experimental 
triple-point temperatures of ^^Ne and ^°Ne is 0.147 K)^ 
while our PI MC simulations predict a value of 0.135 K. 
It is customary to represent the isotopic composition of 
a neon sample of mass m by a variable x, defined as 



where 77120 =20 amu. The PI MC results for the triple- 
point temperature difference with respect to natural neon 
are presented in Fig. [10] as a function of x. The calculated 
TPIE is compared to the experimental data taken from 
Ref. Q- The linear behavior experimentally found as 
a function of x is reasonably reproduced by our simula- 
tions, however our slope is smaller that the experimental 
one. 



V. VAPOR PRESSURE ISOTOPE EFFECT 

The VPIE is experimentally tabulated as the logarith- 
mic ratio of the isotopic vapor pressure, ln(p'/p), where 
p'and p are the vapor pressures of isotopes m and m', 
with m > ml . We have calculated the VPIE for ^^Ne 
and ^°Ne at the critical points of both isotopes. First 
we consider the critical temperature of the lightest iso- 
tope, ^"Ne. Note that this temperature is lower that the 
triple-point temperature of ^^Ne. This implies that at 
this temperature the vapor of ^^Ne may coexist only with 
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Figure 11: Determination of the vapor pressure of ^"Ne and 
^^Ne at the triple point of ^"Ne by the study of the sohd-gas 
coexistence. 



Tjp( Ne) = 24.673 K 
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Figure 12: Determination of the vapor pressure of '^"Ne and 
^^Ne at the triple point of ^^Ne by the study of the liquid-gas 
coexistence. 



Table V: VPIE of ^°Ne and ^^Ne for solid (s) and Uquid (1) 
phases at triple-point conditions, p'and p denote the vapor 
pressure of ^"Ne and '^^Ne, respectively. 



In(p7p)^" 


Hp'/pV 


PI MC 0.056 


0.047 


exper.'^ 0.058 


0.047 



"Values obtained at Ttp{^°Ne) 
''Values obtained at Ttp(22Ne) 
'^Ref. 



its solid phase (as a visual help see the phase diagram in 
Fig. m where there appear two critical points at different 
temperatures). Therefore we have calculated the vapor 
pressure of '^'^Ne and ^ONe at TtpC^^Ne )=24.538 K by 
the study of the solid-gas coexistence. In Fig.lTTjwe have 
plotted the pressure dependence of the Gibbs free energy 
for both solid and gas phases of ^^Ne and ^"^Ne. The free 
energy of the solid was calculated by a RS simulation us- 
ing Eq. ([5]) , while the free energy of the gas was obtained 
from the classical EoS.'^* The vapor pressures of ^°Ne, p' 
= 0.313 bar, and of ^^Ne, p =0.296 bar, are read from 
the crossing points in the figure. The resulting value of 
hi{p' /p) is compared to available experimental data in 
Tab. El 10 

The VPIE at the triple point temperature of the heav- 
ier isotope, ^^Ne, has been calculated by the study of 
the liquid-gas coexistence shown in Fig. [T^ Gibbs free 
energies for the liquid phases were derived from PI MC 
simulations while the gas phase results correspond to the 
Johnson parameterization. '^'^ Here the vapor pressure of 
20Ne and ^^Ne are p' = 0.330 bar and p =0.315 bar, and 
the calculated value of ln{p' /p) shows an excellent agree- 
ment to the experimental data (see Tab. |V|.i"'^^ Our 
simulations correctly reproduce that the VPIE is larger 
in the solid than in the liquid phase. 



VI. CONCLUSIONS 

The calculation of the phase diagram of neon by means 
of path integral simulations demonstrates that free en- 
ergy techniques, previously used in classical simulations, 
can be also effectively employed to study quantum prob- 
lems. In particular, we have made extensive used of the 
Adiabatic Switching2i and Reversible Scaling22 methods, 
that are based on algorithms where either the Hamilto- 
nian, a state variable (pressure, temperature) or even an 
atomic mass are adiabatically changed along a simulation 
run. As long as this change is performed slowly (adiabat- 
ically), the associated reversible work is equal to the free 
energy difference between initial and final states. These 
free energy algorithms can be easily implemented in any 
code prepared for equilibrium simulations. 

Our quantum simulations of the phase diagram of neon 
have focused on the study of two types of effects. The 
first one is a theoretical quantity not accessible by exper- 
imental measurements. This "quantum effect" is defined 
as the difference in phase coexistence properties obtained 
by applying either quantum or classical statistical me- 
chanics. Such effects are not experimentally accessible, 
as there is no way to perform coexistence measurements 
of neon phases in the classical limit. Therefore, the cal- 
culation of these effects is interesting only if we want to 
quantify the errors of a classical treatment of phase co- 
existence. This interest exists because historically the 
approaches used to study phase coexistence of rare gas 
atoms were based on classical statistical mechanics. The 
second studied effect is real, in the sense that it can be 
experimentally measured. We refer here to isotope effects 
in the phase coexistence of neon, such as the triple-point 
temperature or vapor pressure. Classical statistical me- 
chanics incorrectly predicts that phase coexistence is in- 
dependent of the atomic mass, and therefore an isotope 
effect is classically forbidden. 
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We have found that quantum effects in the phase di- 
agram of neon are significant at pressures below 2 x 10'^ 
bar. The whole solid-gas and liquid-gas coexistence lines 
are found below this limit. The main quantum effect 
found in the P — V domain for these two lines is a shift 
of about 1.5 K toward lower temperatures. For the case 
of solid-liquid coexistence, a temperature shift is found 
that decreases with pressure, from a value of 1.5 K at 
triple point conditions to about 0.6 K at 2 x 10'^ bar. No 
significant quantum effects are found in the value of the 
triple-point pressure of neon. 

The employed Lennard- Jones parameters have been 
able to reproduce reasonably the experimental melting 
line of neon in the studied pressure range up to 2 x 10'^ 
bar. However, the calculated solid-gas line shows a 
rigid shift of about 0.6 K toward higher temperatures 
when compared to experimental data. By decreasing the 
value of the Lennard- Jones parameter e it is possible to 
shift the solid-gas line towards lower temperatures, but 
the solid-liquid line is shifted also by the same amount. 
Therefore it seems impossible to reproduce correctly both 
coexistence lines by employing a simple Lennard- Jones 
pair potential. The consideration of three body terms 
in the interaction potential is likely a necessary step to 
achieve better agreement to experiment. This considera- 
tion applies also for the liquid-gas coexistence, specially 
when approaching critical point conditions. The triple 



point is reasonable reproduced by our quantum simula- 
tions, in particular the triple-point temperature and co- 
existing densities, while the predicted triple-point pres- 
sure is too low. 

Our study of isotope effects in the triple-point temper- 
ature and vapor pressure of neon shows good agreement 
to experimental data. These isotope effects depend on 
the free energy difference between neon isotopes, that is 
a function of the kinetic energy and therefore not as de- 
pendent on the details of the interaction potential. The 
calculated difference between the triple-point tempera- 
tures of ^^Ne and ^"Ne was 0.135 K, while the experi- 
mental result is 0.147 Kj^ Our quantum simulations cor- 
rectly reproduce that ^°Ne is more volatile than ^^Ne and 
the calculated vapor-pressure isotope effect in both solid 
and liquid phases at triple point conditions show near 
quantitative agreement to experimental data. 
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